Biophysical  Journal  Volume  95  December  2008  5021-5029 


5021 


Free-Energy  Profiles  of  Membrane  Insertion  of  the  M2  Transmembrane 
Peptide  from  Influenza  A  Virus 


In-Chul  Yeh,*  Mark  A.  Olson/  Michael  S.  Lee,*+:|:  and  Anders  Wallqvist* 

*Biotechnology  High  Performance  Computing  Software  Applications  Institute,  Telemedicine  and  Advanced  Technology  Research  Center,  U.S. 
Army  Medical  Research  and  Materiel  Command,  +  Department  of  Cell  Biology  and  Biochemistry,  U.S.  Army  Medical  Research  Institute  of 
Infectious  Diseases,  Fort  Detrick,  Maryland;  and  ^Computational  and  Information  Sciences  Directorate,  U.S.  Army  Research  Laboratory, 
Aberdeen  Proving  Ground,  Maryland 

ABSTRACT  The  insertion  of  the  M2  transmembrane  peptide  from  influenza  A  virus  into  a  membrane  has  been  studied  with 
molecular-dynamics  simulations.  This  system  is  modeled  by  an  atomically  detailed  peptide  interacting  with  a  continuum 
representation  of  a  membrane  bilayer  in  aqueous  solution.  We  performed  replica-exchange  molecular-dynamics  simulations 
with  umbrella-sampling  techniques  to  characterize  the  probability  distribution  and  conformation  preference  of  the  peptide  in  the 
solution,  at  the  membrane  interface,  and  in  the  membrane.  The  minimum  in  the  calculated  free-energy  surface  of  peptide  insertion 
corresponds  to  a  fully  inserted,  helical  peptide  spanning  the  membrane.  The  free-energy  profile  also  shows  that  there  is  a 
significant  barrier  for  the  peptide  to  enter  into  this  minimum  in  a  nonhelical  conformation.  The  sequence  of  the  peptide  is  such  that 
hydrophilic  amino  acid  residues  at  the  ends  of  the  otherwise  primarily  hydrophobic  peptide  create  a  trapped,  U-shaped  conformation 
with  the  hydrophilic  residues  associated  with  the  aqueous  phase  and  the  hydrophobic  residues  embedded  in  the  membrane.  Analysis 
of  the  free  energy  shows  that  the  barrier  to  insertion  is  largely  enthalpic  in  nature,  whereas  the  membrane-spanning  global  minimum  is 
favored  by  entropy. 


INTRODUCTION 

The  membrane  insertion  of  proteins  or  peptides  plays  im¬ 
portant  roles  in  the  mechanisms  of  viral  infections,  toxin 
actions,  and  antimicrobial  defense  (1-3).  The  proper  under¬ 
standing  of  the  mechanism  of  membrane  insertion  of  peptides 
and  proteins  forms  the  basis  for  developing  therapeutic  in¬ 
terventions  against  bacterial  and  viral  diseases  (4,5).  Since 
typically  the  biologically  active  membrane-protein  or  peptide 
functions  only  in  the  membrane  environment,  it  is  critical  to 
have  an  understanding  of  the  membrane-bound  structure. 
However,  despite  significant  advances  in  experimental  tech¬ 
niques,  only  a  limited  number  of  membrane-protein  struc¬ 
tures  have  been  experimentally  determined  (6).  The  same 
experimental  bottlenecks  exist  for  small,  membrane-bound 
peptides.  In  addition,  these  peptides  may  also  undergo  large 
conformational  changes  as  part  of  the  natural  membrane- 
insertion  process.  The  structure  and  dynamics  of  membrane- 
bound  peptides  have  been  partly  investigated  experimentally 
by  studying  the  insertion  process  of  carefully  designed  small 
synthetic  peptides  (7,8).  For  these  types  of  studies,  atomic- 
detailed  molecular  dynamics  (MD)  simulation  techniques 
represent  a  valuable  complementary  methodology  to  inves¬ 
tigate  membrane-insertion  of  peptides  (7).  Thus,  computer 
simulations  of  membrane  insertion  of  peptides  have  been 
performed  based  on  various  models  of  membranes  and  pro- 
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teins  ranging  from  full  all-atom  to  coarse-grained  models 
with  different  levels  of  complexities  (9-14). 

There  are  two  conflicting  views  on  the  mechanism  of  the 
membrane  insertion  of  peptides.  One  view  is  that  the  peptide 
folds  from  an  unstructured  solvent  state  to  a  helix  at  the 
membrane  interface  before  the  insertion.  This  avoids  the  high 
energetic  cost  of  desolvating  the  hydrogen  bonds  of  the 
peptide  backbone  in  solution  before  the  peptide  associates 
with  the  hydrophobic  membrane  (8).  The  other  view  holds 
that  the  peptide  will  fold  into  a  helix  once  inside  the  mem¬ 
brane.  This  view  is  based  on  the  results  of  a  replica-exchange 
MD  (REMD)  simulation  (15)  of  a  WALP16  synthetic  model 
peptide  with  explicit  solvent  and  lipid-bilayer  molecules, 
where  interface  folding  was  not  observed  (11).  This  mode  of 
insertion  and  intramembrane  folding  was  accompanied  by  a 
large  increase  in  the  system  entropy,  which  compensated  for 
the  desolvation  penalty  (11).  The  membrane-insertion  mech¬ 
anism  of  the  individual  peptide  may  be  dependent  on  its  amino 
acid  sequence  composition.  However,  it  is  not  sufficient  to 
predict  the  preferred  membrane-insertion  mechanism  of  the 
individual  peptide  with  the  sequence  data  alone.  A  more  ex¬ 
tensive  study  with  a  diverse  set  of  model  peptides  is  needed  to 
reach  a  general  conclusion  on  the  mechanism  of  the  membrane- 
insertion  process. 

The  M2  protein  from  influenza  A  vims  is  an  essential 
component  of  the  viral  envelope  and  forms  a  four-helix 
bundle  that  exhibits  a  highly  selective,  pH-regulated,  proton 
ion-channel  activity.  The  influenza  A  vims  enters  the  infected 
cell  by  endocytosis,  and  the  interior  of  the  virion  must  be 
acidified  while  it  is  contained  in  the  endosome  as  a  prereq¬ 
uisite  for  uncoating  (release  of  genetic  material  to  the  cyto- 
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plasm)  (16).  The  proton  channel  formed  by  the  M2  proteins 
provides  this  acidification  function  and  is  a  potential  thera¬ 
peutic  target  (16).  The  M2  transmembrane  peptide  (M2-TMP) 
is  a  truncated,  synthetic  peptide  consisting  of  25  amino  acids 
spanning  the  transmembrane  domain  of  the  original  97-amino- 
acid  M2  protein.  The  sequence  of  the  M2-TMP  is  modeled  as 
SSDPLVVAASIIGILHLILWILDRL  (17).  This  sequence  is 
effectively  hydrophobic  but  contains  hydrophilic  anchoring 
residues  aspartate  (D)  and  arginine  (R)  at  the  ends.  Although 
the  truncated  C-  or  N-terminal  regions  of  the  full-length  M2 
protein  undoubtedly  play  important  roles  in  the  viral  life 
cycle,  an  ion  channel  activity  was  demonstrated  with  the 
tmncated  M2-TMP  (16,17).  Recently,  Stouffer  et  al.  (18) 
determined  the  high-resolution  stmcture  of  the  M2-TMP  in 
the  presence  of  an  amantadine-like  inhibitors  with  x-ray 
crystallography.  In  a  companion  article,  Schnell  and  Chou 
(19)  used  NMR  techniques  to  determine  the  channel  structure 
of  a  slightly  differently  tmncated  version  of  the  M2  protein. 
Both  stmctural  assemblies  share  general  structural  features  of 
the  M2  channel,  such  as  a  four-helix  bundle  state.  However, 
they  also  show  significant  differences  in  structural  details 
related  to  the  amantadine  inhibition  mechanism.  Further 
stmctural  studies  need  to  be  carried  out  to  resolve  these  dif¬ 
ferences.  Additionally,  the  M2-TMP  four-helix  bundle  has 
been  studied  by  explicit  atomistic  MD  simulations  with 
carefully  chosen  initial  starting  configurations  to  address  is¬ 
sues  regarding  structure  and  dynamics  of  the  bundle  in  the 
membrane  (20-22).  These  studies  did  not  address  the  sta¬ 
bility  of  the  four-helix  bundle  in  the  membrane  interface.  In 
contrast,  Bu  et  al.  (23)  used  an  implicit  solvent/membrane 
models  to  determine  the  optimal  aggregation  number  of 
helices  in  the  membrane.  Even  though  the  M2-TMP  tetramer 
is  more  relevant  to  the  biological  function  of  the  proton 
channel,  the  M2-TMP  monomer  in  the  membrane  environ¬ 
ment  has  been  studied  experimentally  (24)  and  computa¬ 
tionally  (12-14,25)  to  understand  its  structure  and  dynamics 
in  the  absence  of  multimeric  interactions.  In  these  computa¬ 
tional  studies,  the  M2-TMP  was  either  preinserted  into  the 
membrane  (12,25)  or  its  backbone  was  fixed  to  the  helical 
configuration  (13,14).  However,  missing  from  these  studies 
is  the  complete  picture  of  the  stmctural  and  energetic  changes 
undergone  by  the  M2-TMP  upon  penetrating  the  membrane 
from  an  aqueous  solution. 

Here,  we  report  an  extensive  MD  simulation  study  of  the 
insertion  of  the  M2-TMP  into  a  membrane.  Different  initial 
conditions,  replica  exchange  (15),  and  biased  sampling 
techniques  with  peptide  position  restraints  (umbrella  sam¬ 
pling)  have  been  used  to  ensure  sampling  of  sufficient  con¬ 
formations.  This  is  done  to  avoid  trapping  the  system  in  local 
minima  and  to  ensure  that  enough  conformations  of  the  M2- 
TMP  are  sampled  at  all  locations  across  the  membrane  in¬ 
terface.  Im  and  Brooks  (9)  contemplated  the  use  of  REMD 
simulations  combined  with  umbrella  sampling  to  obtain 
membrane-insertion  free-energy  peptide  profiles.  Recently, 
MD  simulations  employing  an  extensive  two-dimensional 


(2D)  biased-sampling  technique  with  explicit  solvent/mem¬ 
brane  representation  were  performed  to  study  the  thermo¬ 
dynamic  stability  of  a  charged  arginine  in  a  transmembrane 
helix  (26,27).  Similarly,  the  distribution  of  individual  amino 
acids  in  explicit  solvent  and  lipid  bilayer  was  studied  using 
umbrella  sampling  techniques  (28,29).  We  find  that  it  is  es¬ 
sential  to  include  both  biased  and  REMD  samplings  to  obtain 
a  reliable  description  of  the  membrane-insertion  process  for 
the  case  of  the  M2-TMP.  Indeed,  all  simulations  except  those 
that  employed  the  REMD  simulation  combined  with  the 
umbrella  sampling  displayed  deficiencies  in  the  sampling  of 
the  peptide  conformation  of  the  M2-TMP.  To  gain  an  un¬ 
derstanding  of  the  mechanism  of  membrane  insertion  of  the 
peptide,  we  calculated  the  potential  of  mean  force  (PMF) 
(30),  or  the  Helmholtz  free  energy,  of  the  peptide  insertion. 
We  find  that,  at  the  global  free-energy  minimum,  the  M2- 
TMP  is  a  fully  inserted  helix  with  a  tilt  angle  with  respect  to 
the  membrane  normal  closely  matched  to  the  experimentally 
measured  value.  In  addition,  we  find  that  there  is  a  significant 
free-energy  barrier  for  the  nonhelical  peptide  to  insert  into  a 
membrane  from  the  aqueous  phase  commensurate  with  the 
traditional  view  of  peptide  insertion.  The  constructed  free- 
energy  surface  also  reveals  an  entropic  barrier  for  the  transfer 
from  the  aqueous  solution  to  the  fully  formed  helix  spanning 
the  membrane. 

MATERIALS  AND  METHODS 

MD  simulations  were  performed  with  the  CHARMM  MD  simulation  pro¬ 
gram  version  c31b2  (31).  We  used  the  all-atom  CHARMM  22  force  held 
with  CMAP  modification  for  proteins  (32,33).  The  membrane  and  the  water 
solvent  were  represented  by  an  implicit  membrane/solvent  model  implemented 
in  the  “Generalized  Bom  with  a  simple  Switching”  (GBSW)  module  of 
CHARMM  (25).  The  implicit  membrane  model  has  been  successfully  applied 
to  study  membrane  insertions  of  various  model  peptides  (9).  The  membrane 
thickness  of  25  A,  a  membrane  smoothing  length  of  5  A,  and  the  sur¬ 
face  tension  coefficient  of  0.04  kcal/(mol  A2)  were  used  to  represent  the 
dimyristoylphosphatidylcholine  (DMPC)  membrane  unless  specified  other¬ 
wise  (9,25).  The  DMPC  was  chosen  because  it  is  typically  used  in  experi¬ 
mental  studies  of  the  M2-TMP  system  (34). 

Initial  structures  of  the  M2-TMP  for  simulations  were  taken  from  either  a 
fully  extended  conformation  or  an  ideal  a -helix  structure  as  determined  by 
solid-state  nuclear  magnetic  resonance  (NMR)  techniques  (34).  These 
starting  structures  were  optimized  by  energy  minimization  with  a  distance- 
dependent  dielectric  constant  of  £  =  4 r  that  included  50  steps  of  initial 
steepest-descent  minimization  followed  by  minimization  with  the  implicit 
membrane  model  over  200  steps.  During  the  minimization,  CJCp  atoms  of 
all  residues  were  restrained  to  their  initial  positions  with  a  force  constant  of 
0.5  kcal/mol,  except  for  the  last  100  steps.  The  resulting  minimized  con¬ 
formations  were  used  as  starting  structures  of  the  subsequent  MD  simula¬ 
tions. 

A  time  step  of  2  fs  was  used  in  the  MD  simulations.  Langevin  dynamics 
with  a  friction  coefficient  of  5.0  ps1  was  used  for  the  temperature  control 
(9).  Covalent  bonds  between  the  heavy  atoms  and  hydrogens  were  con¬ 
strained  by  the  SHAKE  algorithm  (35).  Distances  used  for  the  onset  of  a 
switching  function  for  nonbonded  interaction,  the  cutoff  for  nonbonded  in¬ 
teractions,  and  the  cutoff  for  nonbonded  list  generation  were  20, 22,  and  25  A, 
respectively.  Coordinates  were  saved  at  every  1  ps  for  further  analysis.  We 
performed  REMD  simulations  (15)  by  running  the  CHARMM  MD  simula¬ 
tion  program  with  the  Multiscale  Modeling  Tools  for  Structural  Biology  (36). 
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Temperatures  used  for  REMD  simulations  were  32  temperatures  exponen¬ 
tially  spaced  over  the  range  from  300  to  800  K  unless  specified  otherwise.  The 
exchange  of  the  conformations  in  neighboring  temperatures  was  attempted  at 
every  500  time  steps  (1  ps)  according  to  the  Metropolis  criterion. 

In  simulations  with  the  peptide-position  biasing,  the  Z-component  of  the 
center-of-mass  position  of  the  peptide  with  respect  to  the  membrane  center 
was  restrained  by  the  harmonic  potential  with  a  force  constant  of  2  kcal/(mol 
A2)  using  the  GEO  command  in  the  Miscellaneous  Mean  Field  Potential 
module  in  CHARMM.  The  Z-axis  is  defined  as  an  axis  oriented  parallel  to 
the  membrane  normal  with  the  membrane  center  as  its  origin,  so  that  the 
membrane  with  the  thickness  of  25  A  used  in  our  study  is  bound  by  two 
planes  at  Z  =  ±12.5  A.  The  applied  force  constant  was  found  to  give 
overlapping  histograms  of  the  peptide  Z-positions  for  all  neighboring 
Z-position  biased  systems.  The  weighted  histogram  analysis  method 
(WHAM)  (37)  described  by  Gallicchio  et  al.  (38)  was  used  to  remove  biases 
due  to  temperatures  or  peptide  positions  and  obtain  PMF  or  free-energy 
profiles  of  membrane  insertion  of  the  peptide.  The  statistical  uncertainties  in 
the  WHAM  analysis  were  also  estimated  by  the  procedure  outlined  in 
Gallicchio  et  al.  (38).  We  calculated  the  PMF  as  — RTlog(P),  where  P  is  a 
probability  distribution  as  a  function  of  the  root  mean-square  distance 
(RMSD),  the  Z  position  of  the  peptide,  or  the  tilt  angle  at  temperature 
^obtained  from  the  WHAM  analysis,  and  R  is  the  ideal  gas  constant  in  units 
of  kcal/mol/K.  The  RMSD  between  two  structures  is  defined  as  the  square 
root  of  the  minimum  average  square  distance  between  respective  backbone 
atoms  of  the  two  structures  with  respect  to  all  rigid  body  rotations  and 
translations  (39).  The  RMSD  was  calculated  with  respect  to  the  a-helix 
structure  determined  by  solid-state  NMR  techniques  (34)  unless  specified 
otherwise.  The  tilt  angle  was  defined  as  an  angle  between  the  long  axis  of  the 
peptide,  as  determined  from  the  inertia  tensor,  and  the  membrane  normal. 

The  durations  of  the  simulations  performed  in  this  study  were  as  follows: 
The  REMD  simulation  of  the  M2-TMP  in  aqueous  solution  lasted  5  ns. 
REMD  simulations  of  the  fully  extended  M2-TMP  preinserted  in  the 
membrane  and  placed  outside  the  membrane  lasted  30  and  40  ns,  respec¬ 
tively.  We  performed  8 1  5-ns  MD  simulations  of  the  M2-TMP  using  a  helical 
starting  configuration  with  its  center-of-mass  Z-position  restrained  at  loca¬ 
tions  ranging  from  Z  =  0  A  to  40  A  in  0.5-A  increment.  For  the  31  systems 
restrained  between  Z  =  0  A  and  15  A,  we  performed  additional  4-ns  REMD 
simulations.  This  was  done  to  enhance  sampling  of  the  barrier  regions  of  the 
free  energy. 

The  temperature  ranges  used  for  unrestrained  REMD  simulations  of  a 
fully  extended  M2-TMP  in  the  presence  of  the  membrane  were  from  300  to 
500  K  when  it  was  preinserted,  and  from  300  to  800  K  when  it  was  placed 
outside  the  membrane.  The  average  acceptance  ratios  for  the  replica  ex¬ 
change  during  the  last  10  ns  of  production  runs  were  77%  and  57%  for  the 
narrower  and  wider  temperature  ranges,  respectively.  The  numbers  of  round 
trips  of  replicas  from  the  lowest  temperature  to  the  highest  and  back  were  138 
with  the  average  round-trip  time  of  3.4  ns  for  the  30-ns  REMD  run  with  the 
preinserted  M2-TMP,  and  32  with  the  average  round-trip  time  of  9.1  ns  for 
the  40-ns  REMD  run  of  the  M2-TMP  started  outside  the  membrane.  All 
replicas  were  not  distributed  equally  among  the  temperatures.  However, 
every  replica  visited  each  temperature  at  least  23  times  with  the  preinserted 
M2-TMP,  and  all  the  temperatures  except  for  the  lowest  two  were  visited  at 
least  three  times  by  every  replica  when  the  M2-TMP  was  placed  outside  the 
membrane. 


RESULTS 

We  present  our  results  in  the  following  order:  We  first  report 
on  the  solution  structure  of  the  M2-TMP  using  REMD  sim¬ 
ulations.  This  is  followed  by  the  results  of  REMD  simula¬ 
tions  to  characterize  the  M2-TMP  in  the  presence  of  a 
membrane.  We  then  report  on  the  calculation  of  the  free- 
energy  surface  of  the  M2-TMP  across  the  membrane  inter¬ 


face  using  extensive  REMD  simulations  combined  with  bi¬ 
ased  sampling  techniques. 

M2-TMP  in  aqueous  solution 

To  understand  the  structural  preferences  of  the  M2-TMP  in 
aqueous  solution,  we  performed  an  REMD  simulation  of  the 
peptide  in  the  absence  of  a  membrane.  The  water  solvent  was 
described  implicitly  with  the  GBSW  model.  This  simulation 
was  started  from  the  helical  conformation  shown  in  Fig.  1  a 
and  lasted  4  ns.  The  distribution  of  RMSD  of  structures  at  the 
room  temperature  of  300  K  with  respect  to  the  initial  helical 
conformation  is  shown  in  Fig.  1  b.  This  was  calculated  by 
performing  the  WHAM  analysis  on  the  data  from  the  last  2  ns 
of  the  simulation.  Snapshots  and  the  RMSD  distribution  of 
the  M2-TMP  at  room  temperature  shown  in  Fig.  1,  a  and  b , 
indicate  that  the  initial  helical  structure  of  the  M2-TMP  de¬ 
cays  rapidly  and  is  not  present  in  any  measurable  quantity  in 
the  solution  phase.  To  examine  the  structures  of  the  M2-TMP 
in  aqueous  solution  in  more  detail,  we  calculated  the  distri¬ 
bution  of  all  pairwise  RMSD  between  400  stmctures  ran¬ 
domly  selected  from  the  last  2  ns  of  the  300  K  trajectory.  This 
distribution  is  shown  in  Fig.  1  c,  and  its  lack  of  uniformity  is 
indicative  of  the  presence  of  diverse  solution  structures.  In 
particular,  the  distribution  has  an  isolated  peak  at  a  pairwise 
RMSD  of  0.7  A,  indicating  that  there  may  be  a  set  of  dom¬ 
inant  structures.  Indeed,  it  was  found  that  199  of  the  ran¬ 
domly  selected  400  structures  belong  to  a  set  of  stmctures 
within  a  pairwise  RMSD  of  2  A.  This  set  includes  the  last 
configuration  shown  in  Fig.  1  a.  In  this  set  of  structures,  the 
helicity  is  broken  in  the  middle  of  the  peptide  sequence.  The 
structure  of  the  M2-TMP  at  3.625  ns  in  Fig.  1  a  shows  a 
disruption  of  the  helicity  at  a  slightly  different  location  of  the 
peptide  chain  and  belongs  to  the  second-largest  cluster  with 
38  structures  whose  pairwise  RMSD  are  within  2  A.  How¬ 
ever,  there  is  a  significant  population  of  structures  (163  out  of 
400)  that  are  not  part  of  these  two  main  conformational 
groups,  which  indicates  a  diverse  set  of  peptide  conforma¬ 
tions  in  solution.  This  is  consistent  with  the  experimental 
observations  that  many  a -helical  membrane  peptides  can 
exist  in  a  variety  of  semi-  or  unstructured  forms  in  aqueous 
solutions  (40,41). 

M2-TMP  in  the  presence  of  a  membrane 

We  performed  MD  simulations  of  the  M2-TMP  in  the  pres¬ 
ence  of  a  membrane  with  different  initial  conditions.  However, 
we  determined  that  simple  applications  of  MD  simulations  at 
constant  temperatures  are  not  adequate  to  properly  and  effi¬ 
ciently  sample  intramembrane  folding  or  membrane  insertion 
of  the  M2-TMP.  Therefore,  we  performed  two  REMD  simu¬ 
lations  starting  from  the  fully  extended  peptide  conformation: 
one  preinserted  perpendicularly  into  a  membrane,  and  the 
other  placed  parallel  outside  the  membrane  as  shown  in  Fig. 
2,  a  and  b ,  respectively.  In  the  REMD  simulation  of  the 
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t  =  0  ns  3.625  ns  4  ns 


FIGURE  1  (a)  Snapshots  of  the  M2  transmembrane  peptide  (M2-TMP) 

(sequence:  SSDPLVVAASIIGILHLILWILDRL)  from  the  REMD  simula¬ 
tions  in  the  implicit  solvent  at  300  K.  The  initial  [time  (t)  =  0  ns]  helical 
conformation  corresponds  to  the  experimentally  determined  structure.  ( b ) 
The  RMSD  distribution  of  the  M2-TMP  in  aqueous  solution  without  the 
membrane  obtained  from  the  WHAM  analysis  on  the  last  2  ns  of  trajectories 
from  the  REMD  simulation.  The  RMSDs  are  calculated  with  respect  to  the 
initial  helical  conformation,  (c)  The  distribution  of  pairwise  RMSD  of  400 
structures  randomly  selected  from  the  last  2  ns  of  the  300-K  trajectory.  This 
distribution  was  calculated  by  forming  all  possible  pairs  between  the  400 
structures  and  calculating  the  RMSD  for  each  pair.  Almost  half  of  the 
structures  (199  out  of  400)  belonged  to  a  cluster  in  which  all  members  were 
within  an  RMSD  of  2  A  from  each  other.  This  set  includes  the  last 
configuration  of  the  M2-TMP  at  4  ns  shown  in  a.  The  second-largest  cluster 
(38  members)  contains  the  configuration  of  the  M2-TMP  at  3.625  ns,  shown 
in  a. 


preinserted  M2-TMP,  we  employed  32  temperatures  ranging 
from  300  to  500  K  and  a  membrane  smoothing  length  of 
0.6  A.  The  temperature  range  and  the  membrane  smoothing 
length  were  chosen  to  match  those  used  in  a  previous  REMD 
simulation  by  Im  et  al.  (25)  of  the  preinserted  M2-TMP. 
Ulmschneider  et  al.  (12)  performed  a  replica-exchange  Monte 
Carlo  simulation  of  the  same  system,  but  with  a  different 


FIGURE  2  Two  different  starting  configurations  of  the  fully  extended 
M2-TMP.  (a)  The  M2-TMP  is  preinserted  perpendicularly  into  the  membrane. 
(b)  The  M2-TMP  is  initially  placed  outside  but  parallel  to  the  membrane 
surface. 


implicit  membrane  model.  Our  simulation  lasted  30  ns,  and 
the  last  10  ns  of  the  simulation  were  used  for  the  analysis.  It 
is  important  to  note  that  extensive  equilibrations  are  required 
in  REMD  simulations  to  avoid  simulation  artifacts  (42).  The 
PMF  surface  was  constructed  as  a  function  of  tilt  angle  and 
the  RMSD  from  the  helical  structure.  This  is  shown  in  Fig. 
3  a.  Similarly,  the  PMF  surface  was  also  calculated  as  a 
function  of  the  Z-distance  of  peptide  from  the  membrane 
center  and  the  RMSD  from  the  helical  structure,  as  shown  in 
Fig.  3  b.  The  most  probable  tilt  angle  in  the  PMF  surface  in 
Fig.  3  a  is  consistent  with  the  experimentally  observed  value 
of  38°  ±  3°  in  the  DMPC  membrane  (34).  It  corresponds  to 
the  tilt  angle  needed  to  localize  the  charged  groups  of  the 
helix  to  the  membrane  interface.  The  tilt  angle  is  also  in 
agreement  with  the  35°  ±  2°  observed  in  the  recently  de¬ 
termined  crystal  stmcture  of  the  M2-TMP  four-helix  bundle 
(18).  The  PMF  surfaces  in  Fig.  3,  a  and  b ,  and  examination  of 
the  configurations  of  the  M2-TMP  at  the  lowest  temperature 
show  that  REMD  simulations  of  the  preinserted  peptide 
converged  to  those  of  membrane- spanning  helical  confor¬ 
mation  (RMSD  =  0  A),  consistent  with  results  of  previous 
replica-exchange  simulations  (12,25). 

In  an  attempt  to  simulate  spontaneous  membrane  insertion 
of  the  M2-TMP  from  the  aqueous  solution,  we  also  per¬ 
formed  an  REMD  simulation  of  the  M2-TMP  initially  placed 
45  A  away  from  the  center  and  parallel  to  the  membrane  in  an 
initially  fully  extended  peptide  conformation,  as  shown  in 
Fig.  2  b.  In  this  REMD  simulation,  we  employed  32  tem¬ 
peratures  ranging  from  300  to  800  K  and  a  membrane 
smoothing  length  of  5.0  A.  The  same  set  of  parameters  was 
used  in  REMD  simulations  of  membrane  insertions  of  other 
peptides  with  the  same  implicit  membrane  model  (9).  The 
simulation  lasted  40  ns,  and  configurations  from  the  last  10  ns 
were  used  for  the  analysis.  We  calculated  a  PMF  surface  as  a 
function  of  the  Z-position  of  the  peptide  and  RMSD  with 
respect  to  the  helical  conformation  using  the  WHAM  anal- 
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FIGURE  3  PMF  surfaces  at  300  K  calculated  from  the 
conformations  generated  by  REMD  simulations  of  the  M2- 
TMP.  (a)  PMF  surface  as  a  function  of  the  RMSD  and  tilt 
angle  of  the  M2-TMP  initially  preinserted  into  the  mem¬ 
brane.  ( b )  PMF  surface  as  a  function  of  the  RMSD  and 
Z-position  of  M2-TMP  initially  preinserted  into  the  mem¬ 
brane.  A  membrane-spanning  configuration  of  the  peptide 
is  shown  in  atomic  detail  at  right.  The  hydrophilic  residues 
are  drawn  with  atomic  van  der  Waals  spheres  and  are 
shown  located  at  both  sides  of  the  membrane,  (c)  PMF 
surface  as  a  function  of  the  RMSD  and  Z-position  of  the 
M2-TMP  initially  placed  outside  of  the  membrane  at  a 
distance  of  45  A  away  from  the  membrane  center.  The 
minimum  of  the  PMF  is  set  to  0  by  normalizing  the 
probability  distribution  with  the  maximum  probability. 
Here,  a  representative,  trapped  peptide  conformation  is 
shown  with  the  hydrophilic  residues  at  both  ends  of  the 
peptide  at  the  interface  and  the  hydrophobic  middle  resi¬ 
dues  solvated  in  the  membrane  core. 


ysis.  The  PMF  surface  along  with  the  last  configuration  of  the 
M2-TMP  at  the  temperature  of  300  K  is  shown  in  Fig.  3  c.  In 
this  case,  however,  structures  generated  by  the  REMD  sim¬ 
ulation  did  not  contain  stmctures  that  were  fully  inserted  into 
the  membrane  (Z  =  0  A),  as  seen  by  the  lack  of  any  con¬ 
formations  in  the  lower  part  of  the  PMF  surface  in  Fig.  3  c. 
Instead,  the  REMD  simulation  generates  bent  conformations 
where  the  middle  section  of  the  peptide  is  inserted  into 
membrane  while  the  two  end  groups  remain  near  the  interface 
outside  the  membrane  core,  as  shown  by  the  van  der  Waals 
spheres  at  right  in  Fig.  3  c.  This  trapped  conformation  of  the 
M2-TMP  resembles  the  final  solution  stmcture  shown  Fig. 
1  a  with  an  RMSD  of  3.5  A.  Both  structures  share  a  common 


motif  of  a  U-shaped  conformation  even  though  the  details 
of  the  stmctures  differ.  However,  U-shaped  conformations 
of  the  M2-TMP  in  the  membrane  interface  and  the  solution 
are  thought  to  have  different  origins.  The  U-shaped  con¬ 
formation  of  M2-TMP  in  the  solution  may  arise  from  the 
interresidue  hydrophobic  interactions.  The  U-shaped  con¬ 
formation  at  the  membrane  interface  is  mainly  due  to  the 
partitioning  of  hydrophilic  end  groups  and  hydrophobic 
groups  between  aqueous  and  membrane  environments,  re¬ 
spectively.  According  to  hydrophobicity  scales  developed  by 
Wimley  and  White  (43)  and  more  recently  by  Hessa  et  al. 
(44,45),  both  N-  and  C-terminal  ends  of  the  M2-TMP  contain 
strongly  hydrophilic  residues  Asp  and  Arg,  whereas  the 


Biophysical  Journal  95(11)  5021-5029 


5026 


Yeh  et  al. 


middle  section  of  the  M2-TMP  includes  hydrophobic  resi¬ 
dues,  e.g.,  Leu,  Val,  lie,  and  Trp.  Although  these  strong 
hydrophilic  amino  acid  residues  may  play  important  roles  in 
anchoring  this  peptide  in  specific  position  and  tilt  angle, 
they  also  appear  to  present  a  significant  barrier  to  the  full 
membrane  insertion  of  this  peptide  in  the  current  implicit 
membrane  model.  We  also  did  not  observe  a  spontaneous 
peptide  insertion  in  an  REMD  simulation  with  the  different 
smoothing  parameter  of  0.6  A  (results  not  shown).  Our 
observations  on  the  lack  of  spontaneous  insertion  are  thus 
not  strongly  parameter  dependent.  A  remaining  question  is 
whether  we  can  identify  the  global  conformation  minimum 
from  these  two  different  peptide  conformations  shown  in 
Fig.  3.  The  results  up  to  this  point  suggest  that  a  simple 
application  of  REMD  simulation  with  temperature  pertur¬ 
bation  is  not  sufficient  to  overcome  the  problem  of  initial 
condition  dependency  and  determine  the  global  minimum 
structure  of  the  peptide. 

Free-energy  REMD  simulations  with 
peptide-position  biasing 

To  avoid  trapping  the  M2-TMP  in  local  minima  of  the  free- 
energy  surface  and  to  sample  peptide  conformations  at  a 
wider  range  of  points  across  the  membrane/water  interface, 
we  performed  MD  simulations  of  the  M2-TMP  with  its 
center-of-mass  Z-position  restrained  at  multiple  locations 
across  the  membrane  interface.  These  simulations  were 
started  from  the  helical  conformation  based  on  the  NMR 
experimental  structure.  The  peptide  was  initially  oriented 
parallel  to  the  membrane  interface.  The  center-of-mass 
Z-position  of  the  peptide  was  restrained  at  81  points  across 
the  membrane  interface  in  intervals  of  0.5  A  starting  at  the 
center  of  the  membrane  (Z  =  0  A)  and  ending  at  Z  =  40  A. 
Only  the  upper  half  of  Z-positions  (Z  >  0  A)  was  considered 
because  the  peptide-membrane  interaction  is  expected  to  be 
symmetric  with  respect  to  the  membrane  center.  These 
MD  simulations  with  peptide  position  restraints  lasted  5  ns. 
Analysis  of  these  MD  simulations  revealed  a  lack  of  sam¬ 
pling  of  peptide  conformations  characterized  by  the  de¬ 
pendency  on  initial  conditions.  Therefore,  to  ensure  proper 
samplings  of  peptide  conformations  of  the  M2-TMP  across 
the  membrane  interface,  we  performed  additional  4-ns  REMD 
simulations  starting  from  each  of  the  restrained  M2-TMP/ 
membrane  systems  after  5-ns  MD  simulations.  REMD  sim¬ 
ulations  were  not  extended  to  points  past  Z  =  15  A  since 
peptide  conformations  beyond  this  distance  are  not  relevant 
to  the  membrane-insertion  process.  It  is  to  be  noted  that 
conformations  of  the  M2-TMP  in  the  limit  of  large  Z  were 
sampled  by  the  REMD  simulation  of  the  M2-TMP  in  aque¬ 
ous  solution  as  described  above.  The  last  1  ns  of  REMD 
simulation  trajectory  was  used  for  the  analysis.  We  removed 
the  biasing  imposed  by  the  Z-position  restraints  and  calcu¬ 
lated  PMF  surfaces  as  a  function  of  Z-positions  and  RMSD 
with  the  WHAM  analysis.  The  calculated  free-energy  surface 


in  Fig.  4  a  shows  the  relative  free-energy  change  of  each 
conformation  with  respect  to  the  distance  from  the  membrane 
center  and  the  RMSD  from  the  experimental  helical  confor¬ 
mation.  At  near  Z  =  15  A,  stmctures  with  large  RMSD  values 
are  dominant,  in  agreement  with  the  solution  structure  shown 
in  Fig.  1  b.  However,  near  the  interface  at  Z  =  10  A,  small  and 
large  RMSD  values  coexist  with  similar  free  energies.  This 


RMSD  (A) 


b 


FIGURE  4  PMF  at  300  K  of  the  M2-TMP  calculated  from  REMD  simu¬ 
lations  combined  with  the  center-of-mass  Z-positional  restraints.  The  peptide 
center-of-mass  Z-position  is  restrained  at  various  locations  extending  from 
the  center  of  the  membrane  to  the  point  where  Z  =  15  A  at  0.5 -A  intervals. 
The  final  PMF  surface  is  pieced  together  from  restrained  simulations  by  the 
WHAM  analysis.  ( a )  PMF  surface  as  a  function  of  the  M2-TMP’s  Z-position 
and  RMSD  with  respect  to  the  experimental  helical  conformation.  The 
overall  global  minimum  is  seen  as  the  dark  red  area  in  the  lower  left-hand 
corner  and  corresponds  to  a  helical  conformation  spanning  the  membrane 
with  the  appropriate  experimental  tilt  angle.  ( b )  ID  energy  profiles  obtained 
by  averaging  the  2D  probability  distributions  over  all  RSMD  values.  The 
PMF  profile  as  a  function  of  the  Z-position  of  the  M2-TMP  clearly  shows  the 
overall  hydrophobic  attraction  of  the  solvated  peptide  to  the  membrane,  but 
it  also  shows  a  barrier  to  full  penetration  into  the  membrane  core.  The 
averaged  potential  energy  (A  U)  pro  hie  shows  that  the  partially  inserted 
trapped  peptide  is  stabilized  by  the  enthalpic  component  outside  the 
membrane  core,  and  that  these  states  are  more  energetically  favorable  than 
a  fully  inserted  peptide  at  Z  =  0  A.  The  entropic  contribution  (— T AS)  to  the 
free-energy  surface  was  calculated  from  the  difference  between  the  PMF  and 
the  ensemble  average  potential  energy  (PMF  —  At/).  This  graph  clearly 
shows  that  the  trapped  states  outside  the  membrane  core  are  disfavored 
whereas  the  fully  inserted  peptide  (which  is  experimentally  observed)  is 
favored.  The  statistical  uncertainties  in  the  PMF  and  the  average  potential 
energy  profiles  are  twice  the  standard  deviation  estimated  by  the  WHAM 
error  analysis. 
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surface  shows  that  the  global  minimum  conformation  of  the 
peptide  is  that  of  a  fully  inserted  helix,  as  indicated  by  the 
favorable  free-energy  region  at  small  RMSD  near  the  mem¬ 
brane  center.  However,  the  free-energy  surface  also  shows 
that  there  is  a  significant  barrier,  exceeding  90  kcal/mol,  for 
the  M2-TMP  to  fully  insert  into  a  membrane  in  a  nonhelical 
peptide  conformation.  Indeed,  the  surface  indicates  that  con¬ 
formations  starting  outside  the  membrane  are  funneled  into  a 
secondary  minimum  corresponding  to  a  set  of  nonhelical 
conformations  trapped  outside  the  membrane  core  (Z  ~7  A, 
RMSD  ~9  A).  The  PMF  surface  in  Fig.  4  a  also  indicates 
that  there  are  considerable  kinetic  barriers  between  the  two 
conformations,  i.e.,  the  U-shaped-trapped  and  membrane- 
spanning  conformations.  For  example,  a  stmctural  transfor¬ 
mation  of  the  peptide  from  a  nonhelical  conformation  to  a 
near  helical  conformation  at  a  peptide  center-of-mass  posi¬ 
tion  ofZ~7  A  followed  by  membrane  insertion  is  associated 
with  a  free-energy  barrier  on  the  order  of  50  kcal/mol.  Fig.  4  a 
also  indicates  that  there  are  narrow,  lower-energy  paths  be¬ 
tween  the  two  conformations  with  a  free-energy  barrier  on 
the  order  of  10  kcal/mol.  However,  these  paths  require  spe¬ 
cific  peptide  conformations  at  specific  Z  positions,  and  their 
significance  may  need  to  be  checked  with  further  simula¬ 
tions.  Since  there  is  no  clear  “path”  in  this  energy  surface 
between  the  two  regions  in  the  surface,  the  preceding  REMD 
simulations  were  never  able  to  switch  between  the  confor¬ 
mations.  The  probability  of  traversing  the  barrier  in  a  regular 
MD  simulation  without  the  biased  potential  is  thus  negligi¬ 
ble,  and  traditional  non-free-energy  simulations  will  always 
show  a  strong  dependence  on  the  initial  conformations. 

To  look  at  the  overall  changes  in  energies,  we  constmcted 
PMF  profiles  in  only  the  Z-direction  from  the  one-dimen¬ 
sional  (ID)  probability  distribution  of  the  peptide  Z-position 
obtained  by  averaging  the  unbiased  2D  probability  distribu¬ 
tion  for  all  values  of  the  RMSD  coordinate.  The  free-energy 
profile  along  the  Z  axis  with  the  estimated  error  from  the 
WHAM  analysis  in  Fig.  4  b  clearly  identifies  a  global  mini¬ 
mum  at  Z  =  0  A  and  displays  a  shallow  local  minimum  at 
around  Z  =  7  A.  This  local  minimum  is  reflective  of  the 
M2-TMP  trapped  outside  the  membrane  core  in  a  nonhelical 
conformation  where  the  hydrophobic  middle  part  of  the 
peptide  and  the  charged  residues  at  the  ends  are  attracted  to 
the  membrane  core  and  the  water  phase,  respectively.  As  the 
M2-TMP  is  further  pushed  in,  the  repulsion  increases  due  to 
the  loss  of  water  solvation  of  the  charged  residues,  and  an 
effective  free-energy  barrier  separating  the  global  and  local 
minima  is  created.  This  barrier  is  centered  at  around  Z  =  5  A. 
The  magnitude  of  the  barriers  in  the  2D  free-energy  profile  in 
Fig.  4  a  appears  more  pronounced  because  of  the  detailed 
structural  breakdown  provided  in  the  2D  profile.  Once  the 
peptide  overcomes  this  barrier,  it  readily  finds  the  global 
minimum  at  the  membrane  center  where  the  hydrophobic 
middle  and  the  charged  residues  at  each  end  occupy  the 
membrane  core  and  the  hydrophilic  solvent,  respectively. 
The  potential  energy  profile  along  the  Z-axis  in  Fig.  4  b 


shows  that  the  overall  minimum  of  the  potential  energy  is 
located  at  around  Z  =  7  A  instead  of  at  the  membrane  center, 
suggesting  that  the  trapped  configuration  is  energetically 
favorable.  The  fact  that  the  fully  inserted  M2-TMP  is  at  the 
global  minimum  is  due  to  entropic  components  in  the  free 
energy.  This  is  illustrated  in  Fig.  4  b ,  where  we  show  the 
difference  between  the  PMF  and  the  internal  energy,  i.e.,  the 
entropic  contribution  (— TAS)  to  the  free  energy.  It  is  to  be 
noted  that  this  estimate  of  the  entropic  contribution  does 
not  include  the  contribution  from  the  nonpolar  solvation 
energy  in  the  implicit  solvent/membrane  model.  The  U-shaped 
trapped  conformation  of  the  M2-TMP  is  favored  by  enthalpy 
but  disfavored  by  entropy,  whereas  the  membrane  spanning 
helical  conformation  is  actually  disfavored  by  enthalpy  but 
favored  by  entropy.  The  favorable  entropic  contribution  to 
the  free  energy  was  also  observed  in  a  simulation  study  of  the 
WALP16  peptide  using  an  atomically  explicit  membrane 
model  (11). 

DISCUSSION 

In  this  study,  we  performed  extensive  MD  simulations  of  the 
M2-TMP  at  various  positions  across  the  membrane  interface 
to  study  the  membrane-insertion  process  of  an  atomically 
detailed  peptide.  The  aqueous  solution  and  the  bilayer  mem¬ 
brane  itself  were  treated  in  the  continuum  approximation  using 
a  generalized  Bom  model.  To  enhance  samplings  of  peptide 
conformational  states  beyond  those  trapped  in  local  minima, 
we  performed  REMD  simulations  with  and  without  peptide 
positional  restraints.  REMD  simulations  without  peptide  po¬ 
sitional  restraints  were  able  to  find  the  experimental  helical 
conformation  of  the  M2-TMP  when  it  was  preinserted,  but 
failed  to  show  a  spontaneous  membrane  insertion  when  the 
peptide  was  initially  placed  outside  the  membrane.  Only  by 
performing  REMD  simulations  combined  with  peptide  posi¬ 
tional  restraints  were  we  able  to  identify  the  fully  inserted 
M2-TMP  helix  as  the  global  minimum  state  in  the  free-energy 
surface,  which  is  in  agreement  with  the  experimental  results. 
These  findings  underscore  the  importance  and  difficulties  of 
sampling  enough  conformational  states  to  obtain  a  meaning¬ 
ful  and  statistically  significant  free-energy  surface  from  MD 
simulations. 

There  have  been  conflicting  notions  regarding  the  mech¬ 
anism  of  peptide  membrane  insertion.  The  central  question  is 
whether  the  a -helix  formation  precedes  the  peptide  insertion 
into  a  membrane  or  not.  In  the  broadly  accepted  model  (8), 
protein  folding  at  the  interface  always  precedes  membrane 
insertion  of  proteins.  However,  a  recent  REMD  simulation 
study  (11)  of  membrane  insertion  of  the  WALP16  peptide 
using  explicit  solvent/lipids  indicates  that  the  interfacial 
folding  is  not  required  for  bilayer  insertion.  The  current  re¬ 
sult  from  the  M2-TMP  using  an  implicit  membrane  model 
suggests  that  the  peptide  would  have  difficulty  penetrating 
the  membrane  in  a  nonhelical  conformation  due  to  a  large 
free-energy  barrier.  It  is  possible  that  peptides,  such  as  the 
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M2-TMP,  that  have  strongly  hydrophilic  groups  follow  dif¬ 
ferent  membrane-insertion  mechanisms  than  the  WALP16 
peptide.  It  is  also  possible  that  detailed  interactions  of  the 
peptide  with  explicit  lipid  headgroups  and  water  molecules 
may  be  important  in  the  initial  stages  of  the  peptide  insertion 
to  escape  any  trapped  conformational  state  in  local  minima. 
In  recent  simulation  studies  of  a  charged  arginine  side  chain 
in  a  transmembrane  helix  using  an  explicit  atom  representa¬ 
tion  of  the  lipid  membrane  and  waters,  it  was  shown  that  the 
charged  side  chain  is  significantly  stabilized  by  favorable 
interactions  with  water  molecules  and  lipid  headgroups 
(26,27).  It  is  suggested  that  the  depolarization  effects  are 
overestimated  in  implicit  membrane  models  that  do  not  ac¬ 
count  for  the  flexibility  of  the  membrane  interface,  and  hence 
the  effective  barriers  to  hydrophilic  peptide  penetration  are 
overestimated  in  such  models  (26,27,46).  In  the  explicit  atom 
simulations,  the  fluctuations  of  motions  of  lipids  and  water 
can  create  vacancies  that  can  facilitate  the  membrane  inser¬ 
tion  of  the  peptide.  However,  performing  fully  atomistic 
simulations  and  employing  extensive  sampling  techniques 
still  remain  a  significant  challenge.  The  implicit  membrane 
model  in  conjunction  with  good  sampling  techniques  can  still 
provide  insights  into  the  peptide  membrane-insertion  process 
by  capturing  the  gross  physical  features  of  the  system.  Im¬ 
provements  of  the  implicit  membrane  model  with  varying 
local  dielectric  coefficients  across  the  membrane  interface  (47) 
could  be  used  to  further  improve  the  implicit  representation. 
Possible  alternatives  to  either  the  costly  explicit-membrane 
model  or  simplistic  implicit  membrane  representations  are 
various  coarse-grained  membrane  models  (48-52).  Fine-tun¬ 
ing  of  the  implicit  membrane  model  with  the  results  from 
simulations  with  explicit  membrane/water  would  be  highly 
desirable  (28).  The  formation  of  multimeric  complexes  may 
also  play  an  important  role  by  facilitating  self-solvation  of 
hydrophilic  peptide  groups  in  the  hydrophobic  membrane  en¬ 
vironment.  Further  experimental  and  computational  studies 
with  a  diverse  set  of  peptides  may  be  needed  to  arrive  at  any 
definitive  conclusion  on  the  mechanisms  of  membrane  inser¬ 
tion  of  peptides.  The  results  from  this  study  of  the  membrane- 
insertion  mechanism  demonstrate  that  careful  analysis  of  the 
simulations  with  proper  samplings  of  conformational  states  is 
needed  to  obtain  meaningful  information  from  the  simulations. 
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